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Abstract 

The motion of assemblies of point vortices in a periodic parallelogram can be described by 
the complex position Zj {t) whose time derivative is given by the sum of the complex velocities 
induced by other vortices and the solid rotation centered at Zj. A numerical simulation 
up to 100 vortices in a square periodic box is performed with various initial conditions, 
including single and double rows, uniform spacing, checkered pattern, and complete spatial 
randomness. Point process theory in spatial ecology is applied in order to quantify clustering 
of the distribution of vortices. In many cases, clustering of the distribution persists after a 
long time if the initial condition is clustered. In the case of positive and negative vortices 
with the same absolute value of strength, the L function becomes positive for both types of 
vortices. Scattering or recoupling of pairs of vortices by a third vortex is remarkable. 



1 Introduction 



A motion of point vortices (P Vs) is one of the simple and fundamental issues in fluid dynamics f3^ . There 
have been a number of numerical studies of PVs, whose motion is bounded by a circular wall. Recently, 
stable triangular vortex lattices have been observed in rotating Bose-Einstein condensates of Na atoms 
[4] in which the circulation of vortices is quantized with the same strength and the lattice contains over 
100 vortices. Many experimental reports of superfluid vortex lattices increased in these several years 
show triangular patterns [11 [5l [6] , which can be compared with the stability analysis of vortex lattices by 
Tkachenko [71 [8] in 1966. Conformal theory of irrotational flows and complex functions are invoked and 
the role of the Weierstrass zeta function is crucial. 

This paper supplements the Letter [T] by showing the spatial distribution of PVs for various type 
of initial conditions. Tkachenko's work 40 years ago and recent strong interests by physicists on vortex 
lattices in Bose-Einstein condensates stimulate the author to consider the numerical study of motions of 
PVs with periodic boundary conditions. 

In order to quantify clustering of the point distribution, we apply point process theory used in spatial 
ecology. The computed K and L functions judge the distribution to be clustered or not. This method can 
be applied to other systems which include information of the point distribution in two or more dimensions. 

2 Velocity field by periodic point vortices 

The velocity field by a single periodic PV of strength k = 27r is given by [H (H [1] 



where z — x + iy and the overbar denotes the complex conjugate. The corresponding streamfunction is 
denoted by 



The Weierstrass zeta and sigma functions (C and a) contain two complex parameters loi and uj2, which 
are half periods on the complex plane. The second term including — 7r/[4Im((I'iW2)] expresses the solid 
rotation in order to cancel the circulation on sides of the parallelogram induced by the PV. 
The equations of motions for PVs of strength 27r/ii, fii = ±1, can be described by 
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u + iv — iC{z) — iflz, 



(1) 



iP = -Re\na{z) + n\z\^/2. 
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For a single type of PVs with the total number TV, the index i is i = 1, ■ ■ ■ , N . For positive and negative 
PVs with the same strength, the index is i — 1, • • • , iV+, + 1, • • • , N{= + N_). There are three 
known conserved quantities; the Hamiltonian and two components of the linear impulse [I]. 



3 K and L functions 



The K function used in spatial ecology ^ is defined by 



N N 




(4) 
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where N is the total number of points, A = N/ S is the number density, S is the cLFGcL, Xj IS the position 
of the jth points, and 9{r) is the step function. A function for the edge correction on the right hand side 
of which usually appears in data analysis of spatial ecology, is unnecessary in the present periodic 
case. The L function is given by 



If the distribution of points is complete spatial randomness (CSR), then K = irr^ and L = 0. If the points 
are clustered, L > Q. The uniformly spaced (US) distribution leads to i < 0. It is generally difficult to 
distinguish between CSR and US distributions only if we see raw data of the point distribution. However, 
the L function gives judgment on the degree of clustering clearly. 

We note that the value of L{r) depends on r. For examples, the checkered pattern gives a positive 
value of L{r) for 2r ~ Ig^ where Ig is the size of the smallest square, but L{r) < for 4r ~ Ig. 

4 Numerical simulation of a single type of PVs 

The system ([3]) with N up to 100 is simulated numerically by using Mathematica 5.2. The textbook by 
Trott (2006) [llj gives a number of numerical codes for Mathematica programming, including samples of 
point vortices with other boundary conditions. 

We first consider the following four initial conditions: (I) an infinite row that is a discrete model of 
the vortex sheet, (II) PVs located randomly in checkered patterns, (III) PVs located randomly in the 10 
X 10 subsquares, and (IV) CSR in the unit square. The typical initial and final distributions of PVs are 
shown in Figures 1-4. The corresponding L functions can be found in [T]. 



L(r) = VTW7 
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Figure 1: (I) The positions of PVs at t = 0.002 and t = 0.01 denoting the discrete approximation of the 
vortex sheet. 



Initial Distribution of PVs (CP) 
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0.2 0.4 y 0.6 0.8 1 




0.2 0.4 y 0.6 0.8 1 



Figure 2: (II) The positions of PVs at t = and t = 0.1, denoting the checkered pattern. 
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Initial Distribution of PVs (US) 
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Figure 3: (III) The positions of PVs at f = and t = 0.1, denoting the uniformly spaced pattern. 
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Initial Distribution of PVs (CSR) 




Distribution of PVs at t=0.1 (CSR) 



°b 02 0.4 X 0.6 0^8 ^1 




'0 0.2 0.4 X 0.6 0.8 1 

Figure 4: (IV) The positions of PVs at t = and t — 0.1, denoting the CSR pattern. 



In Case (I), we have an analytical expression of L{r) for a single row. The roll-up of PVs starts 
immediately corresponding to the Kelvin-Helmholtz instability. L > for the row and it persists at the 
final time t — tf. The relative errors of the conserved quantities are confirmed to be below 10~^ at the 
final time in all numerical examples. 



5 Numerical simulation of two types of PVs 

Next, we consider the three cases of positive and negative PVs, (V) the Karman vortex street, (VI) 
positive and negative PVs located alternately in checkered segments (16 subsquares), and (VII) the CSR 
distribution. We can define similarly the Kim and Lim functions for {I, ni) = (+, +), (— , — ), and (+, — ), 
where +(— ) denotes the positive (negative) PVs. 

In Case(V), the pairing of positive and negative PVs appears. This corresponds to the pairing insta- 
bility of the double row of the staggered PVs i3^. The typical final time can be estimated as 10 te, where 
te is the average eddy turnover time [T] . 

Clustering of both types of PVs is remarkable again in Case (VII) . The vortex pair moves linearly at 
a velocity proportional to 1/h where h is the distance of two PVs. Then, the pair is scattered by a third 
vortex, or recoupling occurs. This type of scattering on the unbounded plane was solved analytically by 
Aref [T^. We made a GIF animation of motions of PVs. The motion is quite interesting. It looks like a 
motion of gas molecules, which does not exist in the case of the single type of PVs. 

Another method to quantify the distribution of PVs in squares of the same size [5] is called the 
Quadrat method 10 in spatial ecology. 
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Distribution of PVs at t=0.001 (Karman Vortex Street) 



Distribution of PVs at t=0.01 (Karman Vortex Street) 
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Figure 5: (V) The positions of positive and negative PVs at t = 0.001 and t = 0.01, denoting the model 
of Karman vortex street. Black (white) circles denote positive (negative) PVs. 



initial Distribution of Two Kinds of PVs (Cfieckered Pattern) 
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Distribution of Two Kinds of PVs at t=0.1 (CP) 
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Figure 6: (VI) The positions of positive and negative PVs a,tt = 0.001 and t = 0.1, denoting the checkered 
pattern. 
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Figure 7: (Vll) The positions of positive and negative PVs at i = 0.001 and t = 0.1, denoting the CSR 
distribution. 



6 Conclusions and Discussions 

Clustering of PVs is studied numerically based on point process theory and the L functions are computed 
from the final distributions of PVs. Clustering persists up to 10 eddy turnover time. We also did a 
numerical simulation for slightly perturbed square vortex lattices, which can be regarded as a special 
circumstance of Case (III). It showed inactive motions of PVs compared with the other examples shown 
in this paper. It may imply that the instability of square lattices requires more than 100 PVs, or the 
rate of growth is very small. The author is grateful for Professor T. Yamagata for support through his 
research on fluid dynamics over these several years. 
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